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C> Abstract 

<n : 

The evolution of a superthermal relic component of matter is studied 
on the basis of non-equilibrium model of Universe and the Fokker-Planck 
type kinetic equation offered by one of the authors. 
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1 Introduction 

o'. 

O" 1 The diffusion equation describing cosmological evolution of superthermal par- 

' tides under the assumption of interactions scaling recovery in the range of su- 

perhigh energies was studied in paper [1] : 

og i d 2 fog 



where the dimensionless function b(r) is determined as: 

m ' 

6(r) - / g(T,x)xdx; (2) 
o 

It is necessary to solve the equation ([1} with initial and boundary conditions of 
the form: 

j^; g(0,x)=G(x); (3) 

b ! lim G(t, x)x 3 = 0, (4) 

moreover the function G(x) must satisfy the integral conditions: 

G(x)x 2 dx = 1; (5) 



G(x)x 3 dx = 1. (6) 



OO 
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As the detailed researches have shown, the solution of kinetic equation (JTJ in case 
of b(r) — obtained in [I] becomes correct only at greater times of evolution, 
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therefore it becomes problematical to establish parameters relationship of this 
solution with the initial distribution. In this paper we discuss the evolution of 
system at small times. 

It is necessary to note that the exact solution of the kinetic equation in 
diffusive approximation (JlJ, satisfying the normalization relations ([3]), (JBJ), is an 
equilibrium ulatrarelativistic Bolzman distribution with conformal temperature 



27 , 

fo = ye" 31 . (7) 

2 Numerical model of the initial distribution 

Let's consider the initial distribution analogous to Fermi-Dirac distribution, 
given in the form of infinitely differentiable and integrable function: 

G (x) = , A - , (8) 

where A,£,y - parameters of the initial distribution. These parameters must be 
such that normalization relations ©,((6]) fulfill automatically. Thus, we have 
two algebraic relations on three parameters, solving which through parameter 
y we find: 

oo 

t 3 dt 



o(t-y) _|_ 1 

t 2 dt f i 2 dt 



£(J/) = ^ ; A(y) = „ ^ iv) (9) 



»(*-») + i J e^-y^ + 1 



Slow convergence of improper integrals in ^) leads to necessity of their trans- 
formation to following form that is more suitable for numeric calculations: 



Mv) 



tdt 



o(t-y) -|- 1 



o 



oo y y 

dx f xdx f dx f xdx . . 

x + 1 J e x + 1 y J e~ x + 1 J e x + 1 ' ^ ' 

oo oo 

oo oo 

f t 2 dt 2 f dx 



oo oo y y y 

/xdx f x dx 2 f dx 2 f X ^ X f x (\\\ 

e x + 1 + J e x + 1 +V J e- x + 1 V J e x + 1 + J e x + 1 ' ^ ' 
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. t 3 dt of dx n 2 f x & x 





oo oo y 

+3y 



x 3 dx 3 f dx 



e x + 1 J e x + 1 J e-' x + 1 
oo o 



> f xdx n fx dx f x dx 



V i^TT + 3 V^TT-y^TT- (12) 



Using known representations of Riemann £ - functions and Bernoulli numbers 
(e.g., see 0): 

OO 

/ ^r da: = (1 " 2l ~ n)r(n)c(n); (13) 



,'2n — l o2n- 1 i 

F+T*- "T^"** < 14) 



and introducing dimensionless functions: 

y 

S(n,y) = J — J^y^, (n = 0, 1, 2, ...), (15) 
o 

(5(0, y) = ln(l + e«)/2) we get expressions for (fT0|) -(IT2l): 

J 1 (y)=yln2 + ^+y5(0,y)-5(l,y); (16) 

J 2 (y) - y 2 In 2 + y^ + |c(3) + y 2 S(0, y)- 
2yS(l,y) + S(2,y); (17) 

Mv) = y 3 \n2 + y 2 ^ + j,|c(3) + ^ttS- 

y 3 S(0, y) - 3y 2 S(l, y) + 3yS(2, y) - 5(3, y). (18) 
Thus, by means of introduced functions (|10l) - (fT2")) we find expression for b(r): 



b (y) = A(y) 

OO 

My) f xdx My) My) 



xdx 



o${y)x-y _|_ i 



e(y)J e-y + i Ji(y) 



(19) 
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Presentation of © by means of S(n, y) and Riemann £ - functions makes nu- 
meric computations simpler. The results of integration are shown on Figure 
1,2 and we see that £(y) is monotone increasing, A(y) is monotone decreasing 
function. 
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Figure 1: Plot of £(y) function. 
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Figure 2: Plot of A(y) function. 
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As a result, normalized initial distribution function is defined by one arbitrary 
parameter, y: 

G (x,y)= e .f( y > , (20) 

which controls the the degree of non-equilibrium of initial distribution ([5J , Fig- 
ure 3. 
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Figure 3: Normalized initial distribution Go(x,y) according to parameter y: 
top-down: y = 1, y = 3, y = 10. The initial distribution J5f is shown by dotted 

line. 



3 Diffusion equation and integral conservation 
laws 

As known (see, e.g. [5]), the strict consequences from general relativistic kinetic 
theory in case of elastic collisions are integral conservation laws of particles total 
number and their energy. The solution of diffusion equation ([1]), obtained on 
the basis of general relativistic kinetic equations, also must satisfy these laws. 
Therefore it is necessary to verify this fact. 

Multiplying by x 2 the both sides of equation ([I} and integrating received 
expression in parts on all interval of impulsive variable x, we get: 



— — f x 2 G(x,T)dx — 
dr J 



g + 2^) G 



(21) 
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Assuming later on: 



dG 

limx 2 G; lim a; 2 — = 0; (22) 

x^O x^O dx 

dG 

lim x 2 G = 0; lim x 2 — = 0; (23) 

x — ^ oo z— >oc OX 



we obtain from (l2~Tj) : 

a; 2 G(a;,T)da; = Const. (24) 



According to ((5|) the constant at the right side of (|24|) equals to 1. 

Multiplying now by x 3 the both sides of ([TJ and integrating received expres- 
sion in parts on all interval, we get: 



d_ 

dr 



x G(x, r)dx 



OO 

" / ^ (S + 2b ^ G ) dx ( 25 ) 
o 

Considering (|22|) . ([24| . ((4]) and once more integrating received expression in 
parts, we find: 

oo 

— / x 3 G(x, r)dx = ~ x 2 G|!° + 
dr J lu 
o 

oo oo 

+ 2 J xGdx - 2&(t) y x 2 Gdx. (26) 



Taking into account relations <[2]) ,d5j), and also (|22|) . (p4|) . we finally obtain from 

OO 

x 3 G(a;, r)dx = Const. (27) 

According to ^ the constant at the right side of (j2"7| equals to 1. 

It is enough in order to G{x) function's degree of magnitude satisfy the 
following strong inequalities for realization of relations (|2"2")) , (|24l) .(l4")): 

G(x)\ x ^ < -, G^l^oo < (28) 
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We also note that realization of the energy conservation law is provided with 
the presence of term with b(r) coefficient at the right side of diffusion equation 
((T|). If this term is missing, the process of energy transmitting at small times is 
not considered. 



G 



4 Expansion of diffusion equation by smallness 
of r 

This approximation match up the early stages of universe evolution, when par- 
ticle interactions are inessential: 



T < 1. 



(29) 



Then the term at the left side of diffusion equation is main and taking into 
account the initial distribution ([3]) we have: 



— =Q,^G(x,t)=G (x) 



(30) 



Expanding the right side of equation ([1}, substituting expression ([317|) as the 
initial distribution and integrating in time variable we obtain the first correction. 
Sequentially iterating this procedure we get the recurring formula for finding 
higher approximations 



/ 1 9 2 

x z ox 



dGk 
dx 



fc / />oo 

2 y~] I Gi j Gk~ixdx 

i-n V JO 



G?T, 



(31) 



Considering that correction of /c-order is proportional to r fc and integrating pip , 
we get: 



-k+l 



G 



1 d 



k+l 



k + 1 x 2 dx 



dx 



i=o V - 70 



g k -ixdx 



(32) 



Differentiating in expression (I32j) , we finally get: 



k + 1 1 9x 2 x <9a; 



i=0 



dx x 



bk-i 



(33) 



We also find the recurring formula for determining function b(r) in equation 
([T|). Considering (f32j) we obtain according to ([2]): 



f 00 da; a 



fc + 1 7 x dx 



dg k 
dx 



Here and further, if not mentioned especially, k G N. 
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<=0 v Jo 



g k -ixdx 



Integrating in parts an integral in ()34|) and supposing henceforward: 
lim xgi — 0; lim x-j^- = 0; 



lim xc^ = 0; lim x-j^- = 0. 

x— >oo a:— foo Ox 



we get finally: 



r fe+l /-oo 



fc + 1 



(34) 

(35) 
(36) 



(37) 



Now we prove that the corrections to normalized initial distribution jS]) can not 
change the total number density and energy density at every step of iterations. 
Multiplying by x 2 the both sides of (|3"!2)) and integrating on all interval, we 

get: 



x Gk+idx 



r k+l 



k + 1 



dx 



i=o V J ° 



g k -ixdx 



Supposing henceforward 

„2 



lim x z gi = 0; lim x 2 -^- = 0; 

x^fO x^O OX 

lim x 2 gi — 0; lim x 2 -j^- = 0. 

a:— >oo a;— >oo 



(38) 

(39) 
(40) 



we obtain from 



x Gk+idx = 0. 



(41) 



Multiplying by x 3 the both sides of (|3^|) and integrating on all interval, we get: 



x 3 Gk+idx — 



r k+l 



k + l 



dgk 
dx 



k / />oo 

2 5Z [ 9i / 9k- t xdx 
<=o V J ° 
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T k+1 
k + 1 



jr- + 2 ^ ( 9i I gk-ixdx 



dx 



dx. 



(42) 



Considering (|39)) . (|40)) . (|27)l and once more integrating received expression in 
parts, we find from (j42|) : 



y x 3 G k +idx = ^ + ^ < - z 2 ,gfc|^ + 2 J xg k dx- 



OO £ 



gk-ixdx I c?x 



After an obvious simplifications 



(43) 



x Gk+idx = 



k + l \ ~ x 9k\ +2 / xgfcdz- 



-2 / xg k dx J x godx- 





-2 I y gk~ixdx J x 2 gidx 



(44) 



and taking into account expressions ([5]), (|59")l , f4"0)) . (|H|1. we get from (|4*4")) at 
last: 



a; Gfc+ida; = 0. 



(45) 



Thus, we certain that the distribution function iterations of each step don't 
change the total number of particles and energy, this is useful tool for calculation 
correctness. It follows from (|4"5|) that the correction of any order is alternating- 
sign on the interval [0, +oo). Hence, small time r approximation is completely 
equivalent to expansion of exact function G(x,t) to Taylor series by r powers. 



4.1 The first order approximation 

As a first approximation according to recurring formula 
Gi = rgi; 



we have: 



91 



]_d_ x 2 (dGo 
x 1 dx V dx 



2b G Q 



(46) 
(47) 
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Substituting expression (jT7|) into © and integrating in parts, also considering 
that function Go(x,y), defined by expression ([8]), with its derivatives approach 
to at x — > oo faster than any of exponential functions and at x —> has finite 
derivatives, we obtain: 



dx ( 8G 



x dx \ dx 
fdG Q 



2b Q G (x) 



\ dx 

dG 
dx 



+ 2b G (x) 
+ 2b G (x) ) dx 



+ 



-tGo(O) +r2b J G dx 
o 

oo 

A(y) e y 



bi(r) = ~T "] ay \ + r2b I G dx. 
e v + 1 J 
o 

Substituting here the initial distribution (jHJ, we find: 



oc 

/ 



G dx = ^4ln(l + e»), 



a (v) 

u 

therefore 

b 1 (r)=TA(y) 



2b ln(l + e y 

a{y) 



e» + 1 



(48) 



(49) 



Then, substituting the function Go(x,y) into (|47p. we determine an explicit 
form of linear approximation to the initial distribution: 
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Figure 5: Evolution of number density of particles, dn(r, x, y) = x 2 G(t, x, y), 
as a first approximation at y = 6. Solid lines from left to right: 
t = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 

, x _ 2A(y) e^ x -y[2bo(y) - g(y)] + 2b (y) 

x {e i(y) X -y + l) 2 
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A(y)ayMy) + 2b (y)}e^-y 



(eZ(v)x-v + 1)2 



- + 



+ 



[e^v) x -v + l) 3 ' 



(50) 



On the figures 4-11 is shown the evolution of distribution of number density 
and their energy density as a first approximation at the values of y = 1, 3, 6, 10, 
when the significant particles part of the initial distribution lies in the area of 
small energy values. 

As we see from these and also following figures: 

1. In number density and their energy density distributions always appear 2 
maximums and 1 minimum, which shifts to area of higher energy values 
with increasing of y parameter; 

2. Always appear such a moment of time, at which the minimum of distri- 
bution function takes negative value; 

3. After a time the first minimum shifts to the area of lower energies, the 
second minimum - to the area of higher energies. 




Figure 6: Evolution of number density of particles, dn(r, x, y) = x 2 G(t, x, y), 
as a second approximation at y = 1. Solid lines from left to right: 
t = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 
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Figure 7: Evolution of number density of particles, dn(r, x, y) — x 2 G(t, x, y), 
as a second approximation at y — 3. Solid lines from left to right: 
t = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 

4.2 The second order approximation 

As a second approximation according to (|32[) we have: 



Substituting then Gn and gi into (|52p . we find the explicit form of the second 
order approximationo 

2 In case that the explicit form of the second order approximation is an unwieldy we don't 
cite it here 



G 2 



(51) 




(52) 
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Figure 8: Evolution of number density of particles, dn(r, x, y) — x 2 G(t } x, y) 
as a second approximation at y = 6. Solid lines from left to right: 



T = 


= 0; 0, 01; 0, 05; 0, 1, dotted line-T = 0, 2. 
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Figure 9: Evolution of number density of particles, dn(r, x, y) — x 2 G(t, x, y) 
as a second approximation at y = 10. Solid lines from left to right: 
t = 0; 0, 01; 0, 05; 0, 1, dotted line-T = 0, 2. 
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Figure 10: Evolution of the energy distribution of particles, 
de(r, x, y) = x 3 G(t, x, y), as a second approximation at y = 1. Solid lines from 
left to right: t = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 
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Figure 11: Evolution of the energy distribution of particles, 
de(r, x, y) = x 3 G(t, x, y), as a second approximation at y — 3. Solid lines from 
left to right: r = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 
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Figure 12: Evolution of the energy distribution of particles, 
de(r, x, y) = x 3 G(t, x, y), as a second approximation at y — 6. Solid lines from 
left to right: t = 0; 0, 01; 0, 05; 0, 1, dotted line-T = 0, 2. 
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Figure 13: Evolution of the energy distribution of particles, 
de(r, x, y) = x 3 G(t, x, y), as a second approximation at y — 10. Solid lines 
from left to right: r = 0; 0, 01; 0, 05; 0, 1, dotted line-T = 0, 2. 
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Figure 14: Evolution of the energy distribution of particles, 
de(r, x, y) = x 3 G(t, x, y), as a second approximation at y — 3. Solid lines from 
left to right: t = 0; 0, 01; 0, 05; 0, 1, dotted line-r = 0, 2. 

5 Conclusion 

Since particles distribution function is nonnegative by definition, it is clear that 
penetration of distribution minimum to the negative values area is an implica- 
tion of corrections smallness conditions failure in the area of concrete energy 
values. Examined approximations sufficiently describe the global properties of 
superthermal particles distribution. The fact of occurrence of two maximums in 
superthermal particles distribution is very significant. As a point of view offered 
in [3] , [4] superhigh energy particles origin model, the first of these maximums at 
a later time can evolve to equilibrium distribution what certify our calculations, 
the second one can give us high-enegry tail of superthermal relic particles. 
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